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ABSTRACT 
The Matrix Interpretive System devised by Wilson and 

modified by Cantin has been transformed into an interactive 
program operable under a time sharing system. Several new 
commands have been added to extend the usefulness of the code 
and give it the capability to offline read, write and punch 


data, to plot graphs and contour maps and to integrate simple 


polynomial expressions. 
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fe LNERODUCT ION 


A. GENERAL DESCRIPTION 

Many engineering problems in such fields as steady-state 
temperature distribution [Ref. 1], stress in elastic structures 
[Ref. 2], and fluid flows [Ref. 3] can'be formulated as problems 
of numerical linear algebra. Until the advent of modern 
digital computers, all but the simplest of these problems had 
to be reserved to specialists for solutions. Even problems of 
moderate size required a considerable computational effort. 
The computer changed all of this. All of the simple matrix 
operations like multiplication, transposition, solution of eigen 
systems can be coded to solve individual problems. The computer 
program described in this report integrates a number of stand- 
ard subprograms of linear algebra [Refs. 4,5,6,7] into a 
problem oriented computer language. This program, heretofore 
called: Matrix Interpretive System was coded in FORTRAN IV. 
The resulting matrix operations sub-language has a syntax and 
rules of its own as well as a certain number of diagnostic 
messages that can be automatically returned to the user of the 
language. The basic specification for the design of the 
miguage structure was that it had to be simple enough to be 
used by persons unfamiliar with computer programming. As 
written, the language could be used to solve stacks of problems 
sequentially; if an error is found in one of the problems a 
message is printed and a solution for the next problem is 
attempted. To be effective in a conversational mode at a 


remote station, the program had to be heavily modified so as 








to leave effective control of the flow of operations to the 
user. This version of the code is called: Interactive Matrix 
Interpretive System. 

The Interactive Matrix Interpetive System is a general 
program written in FORTRAN IV [Ref. 8], it is a problem oriented 
sub-language. The entire program is coded in double precision 
arithmetic. The only subroutines which are used but are not 
contained in the program are a square root and a natural log 


Subroutine. The program uses 160,000 bytes of core storage. 


I. HISTORICAL BACKGROUND 

The original program was developed by Professor E, L. Wilson 
of the University of California in April 1963 for an IBM 7094. 
It was later adapted by Professor G. Cantin of the Naval 
Postgraduate School for a CDC 1604 in February 1966 and for 


the IBM 360/67 in August 1968. 


OF OBJECTIVES 

The objectives of the author's work have been to: 

i, transform the program into an interactive system 
operable under time sharing from a remote typewriter station; 

2. add operations to extend the capabilities of the 
program; 

oc Creace ay Usb 2 MANUAL Stoumacliicate the Use or 


the program. 








II. GENERAL PROGRAM OPERATION 

The Interactive Matrix Interpretive System will be referred 
to as IMIS in this and following sections of this report. IMIS 
1S written in FORTRAN IV language for the IBM 360/67 time 
Sharing system available at the Naval Postgraduate School 
(CP/CMS). A different version called DSMIS is used for batch 
processing. All operations are performed in double precision 
arithmetic. Matrices are stored internally in a single 
dimensional array, A, that is 10,000 double words in sizZe. 
Matrix names are stored in an array called NAME. Arrays, NROW 
and NCOL, are used to store the number of rows and columns of 
each matrix. 

The following is a brief description of the function and 
operation of most of the subroutines of IMIS. Some of the 
subroutines are straightforward and will not be discussed. A 


flow diagram of the program 1S presented in Figure l. 


A. MAIN PROGRAM AND SERVICE SUBROUTINES 
Hh Main Program 
The Main program is used only to initialize internal 


variables and to call the main service subroutine,. SMIS1. 


Z. Main Working Subroutine (SMIS1) 


SMIS1 is the general working subroutine. Some of the 
operations such as duplication of a matrix and envelope of 
a matrix are carried out in SMIS1l. The other operations are 
Carried out in subroutines which are called from SMIS1. Most 


weenie Crror messages are called from SMISI. 








EPRVECE 


Figure l. 
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3. Dis GPngmor va Matiel (hinted) 


MLIST stores each matrix by name, number of rows and 
columns and starting position in the array, A. It checks the 
name of a new matrix to be stored against the names of pre- 
viously stored matrices to see if pyeenene has been used before. 
If there is another matrix with the same name, the previously 
stored matrix will be deleted before the new matrix is stored. 

a. Printing Of EEron Messages (ERRER) 

When an error iS made in an operational command, ERRER 
is called. It prints out an error mesSage, deletes the current 
operational command and returns control to the typewriter. 

5. Locating a Matrix (MFIND) 

Locates the matrix or matrices to be used in an 
operation. It finds the number of rows and columns of a matrix 
and the starting position in the array A. 

o. Calculatingebicgenvelucs and Bigenvectors Nghe yAR) 

Determines all of the eigenvalues and eigenvectors of 
a matrix. It uses the threshold cyclic Jacobi rotation method 
[Ref. 4]. 

7. Formation of Matrix N (P) 
The bicubic polynomial, 
PaC XP¥F4C xP¥A4C XP¥ #C X94C XPV 94C XP XAY +0 xP4 
Cc xy3+c_ xy2+c__ xy+c_ xtc. _y74+c. y24+c _y+c 
9 10 ad fle aes ihe. ns Hes 
is used in the polynomial fit operation. This polynomial 


man be written as 


Jae 








a ees 


where (F,) = ( x2y3 x3y4 x3y,x3 x23 x2? x2y, x4, 
xy3,xy2,xv,x,v3,y¥2,y¥,l ? 


and tC; } = Concer rion a In the POLFIT operation, values 


of P at all of the sixteen nodal points of the grid shown in 
Figure 2 must be used to solve for all of the unknown co- 
efficients of the polynomial using the following matrix 
equation: 


fees =  (N} {Ga} 
: Ce Os, Seo 


where 


> 


(p{l) ) 
uli 


(pl2)) 
(N] : 


Il 


( p(16)) 
1 


Subroutine P forms matrix N for use in the POLFIT operation. 
The bicubic polynomial, P, will be referred to as CUPOL in the 


Following sections. 
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Pigure 2. Grid For Subroutine P 
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Bi OPERATION SUBROUTINES 
1 PNVeT bLNGae wos amie big meme eye idem Sear Ni) 

SYMINV inverts a symmetric matrix, A, by Gauss 
efimination [Ref. 4]. If the matrix is not symmetric, it takes 
the average of the two numbers, A(i,j),A(j,i), and replaces 
both by the average value and prints a warning that the matrix 
was not symmetric. 

2. Adding, Removing and Storing a Submatrix (STOSM) 

STOSM is called to do three operations. 

a) Store a submatrix at an indicated position ina 
matrix when called by the command, STOSM. 

b) The command, RMVSM, calls STOSM to remove a sub- 
matrix from the original matrix. 

c) ADDSM calls STOSM to add a submatrix to a matrix 
at an indicated position. 

Sy, Function Generation (FUNGN) 
This subroutine takes a function described by a set 


of points (X;,Y;) , not necessarily equally spaced, and 
ES ie poe 


calculates, by straight line interpolation, ordinates at equally 


spaced values Gn 


4, Brewed 1 sO cisco S teens) 
PLOT prints a 10 inch by 10 inch graph of 1 to 4 


curves. It automatically scales the graph to the maximum and 


minimum values of X and Y. The axes are plotted as 
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with the characters, X and Y, appearing at every inch on their 
respective axes. When all of the points are plotted on one side 
of either axis, an offset axis replaces the real axis and is 
indicated by 

0 eel A offset X axis 

b) B------ B offset Y axis 
The origin of the axes and the X and Y scale factors are cal- 
culated and printed below the graph. 


2". Determining the Coefficients of a Bicubic Polynomial 
POLFIT) 


The 16 unknown coefficients of a bicubic polynomial 
surface can be determined in terms of 16 known values of the 
polynomial at 16 known values of the independent variables. 


POLFIT solves the matrix equation, 


iN] {c} = [P] 


for the matrix of unknown coefficients, C, where matrix N is 
formed by the subroutine P. Matrix P is a matrix of the values 
of CUPOL at each nodal point and is the input to POLFIT. 

oF Printing Contour Maps (CONTUR) 

CONTUR prints a contour map of CUPOL. There are 21 
possible contour lines which are represented by the alphabetic 
Metters, A through J, 0, Q through Z. The input is either the 
coefficients of CUPOL with the upper and lower limits of X and 
Y or the values of the function at the grid points used for 
POLFIT. In the latter case, POLFIT is called and the coeffi- 


Srents Of CUPOL are calculated. 


eS 








CUPOL is then scanned over the ranges of X and Y to 
determine the maximum and minimum Z ordinate. The range of Z 
is divided into 42 equal parts represented by alphabetic 
characters and blanks assigned to alternating contour values. 
CUPOL is then evaluated at each point on a line, the correct 
character inserted and the line printed. 


(= Homie t HOnsOnldez Oyuc Selinrness Matrix (NORMA) 
FORMK forms the stiffness matrix [Ref. 9], 


| = 


for a simply supported beam having length, L, Young's Modulus, 


E, moment of inertia, I, shear modulus, G, and effective shear 





area, A. 
The stiffness coefficients are: 
_ 2ETI 2+8 7 
Ky L ( #55) = 1895 
eee 5 - 
Ki2 = 251 ( 336 | S| 
where 
6EI 
B = = 
LAG 
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The matrix equation, 


relates end moments, M): Mo, to their corresponding rotations, 


O14 8. as shown in Figure 3. 





Figure 3. Simple Beam 


S Numerical Integration (ITGRAL) 


The subroutine, ITGRAL, is capable of performing 
numerical integration of one and two dimensional functions 
Memcribed by discrete sets of points, either by Simpson's rule, 


or, if the function to be integrated is known, by Gaussian quad- 


rature. 
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9); Numerical sIntegqrationeOpearoet) ©) lt ke hele al 
Equations (RESPON) 


In the study of vibrations of systems with many 
degrees of freedom, the following system of ordinary differential 
equations governs the problem: 
[M] {X} + [C] {x} + [k} {x} = {7} (1) 
where 
[M] is the mass Matrix, 


[C} is the damping matrix 
and 


[K] is the stiffness matrix. 
In general, if 

{X} = [v] {x} (3) 
where matrix [V] is the matrix of undamped mode Shapes and {x} 
is the vector of the time dependent amplitude function for each 
mode, we get after substitution of (2) Si eons: 

[M} [Vv] {x} + [c] [Vv] {x} + [K] [v] {x} = {p} (3) 
After premultiplication by [v]? and the use of the orthogonality 
properties of normal modes: 


fae (M) [Vl = [-m-} 


[~k~] = [~w2m_] 


fees] 


fat (kK) [v) 


Iv]? [c} [Vv] [—2Awm_] 4) 





(8) such an orthogonality condition imposes some 
restrictions on the form of the damping matrix [C]. 
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where 


Ci 


= 2M; Wz 





and using the notation 


[v]? {Pp} = {p} 


we get 
[~m.] {%} + [—2Awm—] {x} + F-w2m] {x} = {p} (4) 
or 
{i} + [-2r\w_] {x} + Few2.] {x} = {p*} (5) 


where p* = Ines 
De P./ ; 


This last set of equations is completely decoupled 
and can easily be numerically integrated. The subroutine 


RESPON performs this task. 
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IIL. OPERATIONAL COMMANDS 

The following is a list of the operations that IMIS con- 
tains. The operations must conform to the standard format 
shown in Figure 4. OPERATION stands for any of the symbolic 
command names mentioned below. A, B, C and D are matrices used 
as arguments of the operation. Nl, N2, N3 and N4 are four 
integer arguments that must always be right justified. S1 and 
S2 are two floating point arguments. The arguments can be 
input or output values and are not all necessary for each 
command. 

There are two versions of the matrix interpretive system 
in use at the Naval Postgraduate School. One is for time 
sharing and one for batch processing. The commands denoted by 
fa 6)6CUCan 0Oonly be used in conjunction with the time sharing 
option. The READ, WRITE and PUNCH commands make use of 
sequential input/output files which are available with the 
FORTRAN IV language [Ref. g}. (b) 

The READ command uses the file FTO4F001. This file is 
used to load matrices into IMIS from punched cards and is read 
into the computer prior to executing IMIS. Any data stored in 
the file is automatically erased when new data is read into the 
file. 

File FTO8F001 is used to sees BAY output which the user 


wants to print offline. This file is not automatically erased 


(b)A certain knowledge of the CP/CMS time sharing system 
is needed when using a remote terminal of the computer. It is 
presumed that the reader is already familiar with this system 
[Ref. 10]. 
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when new data is read into the file. New data iS appended to 
the data that is already stored in the file. Therefore, the 
file must be erased before executing IMIS. To erase the file, 
the command, 

ERASE) FILE, FT08F001 

is used. After exiting from IMIS, the command, 


OFFLINE, PRINTCC, FILE, FTO8FOO1 


will cause the file to be printed on the offline printer. 

The command, PUNCH, uses file FTO7F001 to store data to 
be punched onto cards. This file appends new data at the end 
of old data so the file must be erased before using IMIS. The 
command, 

ERASE 


Bie, ProyPooL 


b b 
erases the file. To offline punch the data stored in the file, 
the command, 


OFFLINE}PUNCH)FILE}FT07F001 


is issued after exiting from IMIS. Examples of the use of 
these files are given in the sample problems. The operational 
commands are: 

START This command causes the elimination of 
all matrices which are in core storage 
and must be the first command used in 
the SolmeloneOn-a problem ee! | sane ciliorts 
is found during execution of a command, 


an error mesSage is printed on the 


Pe 








LOAD 


typewriter, the current command 


eliminated and control returned to the 


typewriter. If an error is found 


during the execution of a problem with 


PSMLS, the current. problem 1Smeauto> 


matically deleted and the remaining 


cards scanned until a new START or 


STOP command is found. 


Matrix Load: 


A 


N1 


N2 


N3 


Name of matrix to be loaded. 
Number of rows of matrix A. 
Number of columns of matrix A. 
Format control indicate: forntne 
daiwa 

The elements of the matrix are 

typed or punched on cards row-wise 

after the LOAD command according 
to one of the options determined 
by the value of N3: 

if N3= 1 or blank, the format is 
(12F6.0). A maximum of twelve 

numbers per line, uSing as many 

lines aS required to write the 
entire matrix must be used. If 

the decimal point is omitted, it is 

assumed to be at the extreme right 


of the field. 


Zo 








PRINT 


SLOP 


ZERO 


De) 2heNSs— 27 etorl2.0)) 6 numbers peu 
line in fields of 12 columns. 

Cy di N3=.3, (4F18 2.0)),9 4) numbers seer 
line in fields of 24 columns. 

d. if N3= 4, (3F24.0), 3 numbers per 
line in fields of 24 columns. 

Matrix Print: 

A = Name of matrix to be printed. 

Nl= Number of lines, which follow the 
PRINT command, to be printed as a 
label for the matrix. Column 13 
of the command card is used for 
Carrelage Control inf DoMisS sonine 
A one in column 13 will cause the 
matrix to be printed on a separate 
OuEPUL page= 

This operation causes the number of 

lines designated by Nl, which follow 

the operation, to be printed as 
remarks. 

This operation ends a run of problems 

and is the last command to be issued 

in a stack of problems. It can be 
used only once per mune 

Formation of a Null Matrix: 


A = Name of null matrix. 
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ADD 


SUB 


MULT 


TRANS 


INVERT 


SYMINV 


STOSM 


Nl= Number of rows of matrix A. 


N2= Number of columns of matrix A. 


Matrix Addition - 


elneee — eo 


Matrix A is replaced by the sum of 


matrices A and B. 


Matrix Subtraction - [A] - [B]) = [A] 


Matrix A is replaced matrix A minus 


Maver B. 


Matrix Multiplication - [A] «+ [B] = [Cc] 


The product of matrix A times matrix 


B is generated and defined by matrix C. 


Matrix Transpose - [B] = ue 


The transpose of matrix A is generated 


and defined by matrix B. 


Matrix Inversion - 


[A] = [a] ~t 


Matrix A is replaced by the inverse 


Of 1tSeueE 


Symmetric Matrix Inversion - [A] = [A] 


Symmetric matrix A is replaced by the 


inverse of itself. 


Seone sclibmatn: x inv -arMatrax: 


Name of large matrix. 


Name of submatrix to be stored. 


Pow number in 
first element 
Column number 


first element 


Zo 


large matrix of 
of submatrix. 
in large matrix of 


OresubmacrLx. 


1 





STODG 


SCALE 


ADDSM 


RMVSM 


Elements of A are replaced by 
elements of B in designated area. 
Store Row Matrix on Diagonal of Matrix: 


A = Name of Square matrix. 


B Name of row matrix to be stored on 
diagonal of matrix A. 

Diagonal elements of A are replaced by 
elements of B. 


SeqlaneMmltaplrcation —, [Al = Sie 


Add Submatrix to a Large Matrix: 


A Name of large matrix. 

B = Name of submatrix to be added 
into large matrix. 

Nl= Row number in large matrix of 
first element of submatrix. 

N2= Column number in large matrix of 
first element of submatrix. 


Extracts Submatrix [B] from [A]: 


A 


Name of large matrix. 

B = Name of submatrix defined as matrix. 

Nl= Row number in large matrix of 
first element of submatrix. 

N2= Column number in large matrix of 
first element of submatrix. 

N3= Number of rows in submatrix. 

N4= Number of columns in submatrix. 

Matrix A is not modified. 
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LOG 


MS CALE 


RMVDG 


EIGEN 


. 


The Log of each Element 

Each element in the matrix A is 
replaced by the natural log of the 
element. 

Scalar Multiplication 

Each element in the matrix A is 
replaced by B times the element where 
B has previously been defined as a 
1 by 1 matrix. 

Extracts Row Matrix B from Diagonal 
of Square matrix A. 

A = Name of Square matrix. 


B 


Name of row matrix composed of the 
diaconal elements of matrix A. 
Matrix A is not modified. 

Eigenvalues and Eigenvectors: 

The eigenvalues and eigenvectors of 
the system, [A] [X]=A[B][X], are 
determined where [A] is a symmetrical 
matrix and [B] is a diagonal matrix 
of positive elements stored as a row 
matrix. 


A Name of matrix A. 


B Name of matrix B. 


C = Name of matrix of eigenvectors 


stored row-wise. 
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SQREL 


INVEL 


DELETE 


Name of row matrix of eigenvalues, 
r. 

The number of eigenvectors to be 
calculated. The ordering of 
eigenvectors and eigenvalues is 
determined by the sign of Nl as 
follows: 

If Nl is positive, eigenvalues are 
arranged in descending order. 

If Nl is negative, eigenvalues are 
arranged in ascending order. 
Eigenvector normalizing. 

If N2 is blank: 

Eigenvectors are normalized so 
that [C][B]{C]) = [1] 

If N2 is 1: 

Eigenvalues are normalized so 


that [A] = [c]? M.-] [Cc] 


Square Root of Each Element 


Each element in the matrix A is 


replaced by the square root of itself. 


Inversion of Each Element 


Each element in the matrix A 1s 


replaced by the reciprocal of itself. 


Matrix Deletion 


The matrix A is eliminated from core 


SeOrage- 
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DUPL Matrix Duplication: 


A = Name of matrix to be duplicated. 


B = Name of matrix defined to be 
TAICHhe heal wi themat nice ae 
ENVEL Envelope Value of Matrix 


The maximum absolute value in each 
row of matrix A is printed along with 
its column number times an interval. 
A = Name of matrix A. 
Sl= Interval between columns. 

RESPON Numerical Integration of a Set of 


Differential Equations: 


{8}42A [or] (x3 + Rw?) {x}={s}P (t) 
with initial conditions, x(0)=x(0)=0. 
A = Name of a row matrix containing 

the diagonal elements of matrix 
[~wiJ. The dimension of A must 
be (lxn) where n is the number of 
equations in the set. 

B = Name of a column Fane containing 
the elements of {S}. The dimension 
of B must be (nxl). 

C = Name of a row matrix containing 
the values of the function P(t) 
evaluated at equal intervals, At: 


PO) Piao eure «er (mane) 


Z9 








* $2= 


The dimension of C depends on the 
choice of At which must be small 
for accurate results. The size of 
At should be determined by actual 
numerical experimentation with 

the command. 

Name of the output matrix contain- 
ing the calculated values of the 
variables as shown in Figure 5. 
Output interval 1. The values of 
the n unknowns, {x}, are printed 
at iAt,2iAt,3iAt,.-.: af i=], 

all of the values are printed. 
Value of the damping ratio, i. 
Time increment, At, used in the 


numerical integration. 


Xq (1At) , x, (21At) , x, (31At), ree ae, 


Xo (1At) ,xX5 (21At) , x, (31At), ery ae 


x, (LAt) ,x, (21At) ,x, (31At), « es 


Fig. 5. OUTPUT MATRIX, D, FOR RESPON 
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FUNGN 


FORMK 


Function Generation: 
A = Name of matrix which defines a 


function in terms of line segments 


B = Name of matrix to be formed which 
is composed of the y ordinates at 
equal x ordinate intervals of the 
function defined by matrix A. 

Nl= Number of y ordinates to be 
generated. The first ordinate 
will be equal to Ya: 

Sl= X ordinate interval. 

Forms a 2x2 Element Stiffness Matrix: 

A = Name of matrix to be formed. 

The line following the command is as 

follows: 

a. Columns 1-12 Moment of inertia of 
member, I. 

b. Columns 13-24 Effective shear 
area of member, A. 

¢. Columns 25-36 Length of member, L. 

ad. Columns 37-48 Modulus of 
elasticity of member, E. 

The shear modulus, G, is assumed to be 


equal to E/2.4. 


oul 








* 


x 


a 


PUNCH 


READ 


CONT 


WRITE 


Punches Matrix onto Cards: 

A = Name of matrix to be punched. 

Nl= Number of comment cards to be 
punched as a heading for matrix A. 

After the PUNCH command is issued, the 

next line contains the format elected 

by the user as follows: 

(NNTRR.SS) 

where 


NN= Number of elements per card. 


j 


Field width per element. 

SS= Number of Significant digits. 

T = Type of data field. T must be a 
standard F or D Fortran EO EMeree. 

Loads Matrices from file FTO4F001 

The load operation card must have a 

l. in Sl. The last card must be a 

Continue Operation card. 

Continue Operation 

Sl= 1. 

Indicates that all matrices in file 

FTO4F001 have been loaded into IMIS. 

Writes Matrix into File FTO8FOOI1: 

A = Name of matrix to be written. 


Nl= Number of lines of remarks to be 
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PLOT 


wWricteneas a label for thesnatr ue 


A. 
N2 must be equal to l. 
Plots One to Four Curves Described by 
Sets of Points Contained in Matrices 
Pia Cel: 
A = First matrix to be plotted on one 
graph. Plotted as a "+", 
B = Second matrix to be plotted on 
One rgracn. Piletved asa "ee, 
C = Third matrix to be plotted on one 
Gaaph. Plotted as a, 'O 2 
D = Fourth matrix to be plotted on one 
graph. Plotted asa ".". 
SSL, begin jeyete ce Norvciieteepe, (Gries 
will be stored in file FTO8FO001 
N4= Number of lines of remarks to be 


printed at bottom of graph. 


Matrices to be plotted do not have to 


have 
have 
with 
ates 
ordi 
one 


may 


the same number of points. They 
to be entered as N by 2 matrices 
the first column the X co-ordin- 
and "ENC m=Ssecone co UM ee ms anGOn 
nates. When plotting more than 
curve on the same graph, a conflict 


arise if more than one curve passes 
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POLE £2. 


through the same point. The symbol 
printed at that point is determined by 
Ene iellowinag order: Doe.) Beancun. 
Determines the Coefficients of CUPOL: 


A = Input matrix. (18x11) 


B = Output matrix of coefficients of 
CULOL, (Cy 1Cor-+-s1Cy) 
where 
Pe) ee 
2) =H 
A(2) y 
4 —- A4 —— A5-———-A6— — }X 
Be 
FY is —— ls 
Ale —— ans iia 
i Pe NG seal a ee 
Y 
and 


38 ou 3 3 
CUPOL=C)X Y +CoX Y +C3X Y+C4X + 
28 2y 2 2 2 
OE eienG eae Cee 


3 2 
CoXY +Cy AY +Cy  AYHC, Xt 


¥2+C__Y+C 


3 
rerCy 4 15**16 


C13 


34 








CONTUR Contour Maps CUPOn: 
N2= Number of lines of comments. 
N3= Input type indicator. 
a. N3= Blank or l. 
A = oe matrix. Same as for 
POLFIT command. (18x11) 
B= Output Matrix Of Ccoeciirer1enes 
of CUPOL. (16xl) 
Db. (Ne —e2 
Aj= Anput macrix, “(20x1) 
A(1l) = Lower bound of X. 


A(2) 


Upper bound of X. 
A(3) = Lower bound of Y. 


A(4) 


Upper bound of Y. 
A(5) through A(20) are the 
coefficients of CUPOL. 

* N4= Offline print indicator. Stores 
contour map in file FTO9FOO1. 

L[TGRAL Integration of Functions Represented 

by Discrete Values or Polynomials: 

A = Input matrix. (Rxl) 

N3= Integration type indicator: 

a. If N3= 1, Simpson's Rule for one 
dimension. R must be even. 
A(1) = Hy 


A(2) through A(R) are evenly 


Spaced values of function. 
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If N3 = 2, Simpson's Rule for two 


dimensions. 
A(1) = Hy 
202) Hy 


A(3) through A(R) are values of 
function entered row-wise as for 
POEP AE. 

R = N“+2, N must be odd. 
There must be the same number of 


values in both directions. H. and 


x 
Hy Can Cerume iia ie. 

Lf N3= 3, Gauss Quadrature for one 
dimension. 

A(1) = Lower bound of X. 

A(2) = Upper bound of X. 

A(3) through A(14) are the co- 


efficients of an llth degree 


polynomial: 


Vl 10 g 8 7 6 
A3X”"+A,)X" +A, X° +A X +A5X' +AQX + 


> 4 3 2 
AgX +A) 9% +A, 7% +A] 9% +A) 3XtAy 4 


If N3 = 4, Gauss Quadrature for two 


dimensions. 
A(1l) = Lower bound of X. 
A(2) = Upper bound of X. 


A(3) Lower bound Ot. 
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A(4) = Upper bound of Y. 
A(5) through A(20) are coefficients 


of CUPOL. 
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IV. LIMITATIONS 


Ds SIZE AND NUMBER OF MATRICES 

The storage available for all of the matrices in any 
Problem is 80,000 bytes or 10,000 double words. This limit 
1s arbitrary and controlled only by the total core storage 
available at any installation. For the naval Postgraduate 
School installation, 10,000 double words was selected in order 
to reduce turn around time in batch processing and reduce the 
Size of the data set required for time sharing. The total 


mmmeer of matrices which can be stored is 100. 


I. MATRIX NAMES 

Each matrix is referred to by its name composed of one to 
Six alphameric characters including blanks. The name is select- 
ed by the user when the matrix is first loaded or generated by 
an operation. If the name for a previously defined matrix is 
used later as a name of a new matrix, the former matrix is 
deleted before the operation is executed and is replaced by the 
new matrix. This is done so that there is only one matrix with 


a certain name at any one time. 


i. ELOT 

There are 200 points in the Y direction on the graph. The 
X direction has 130 points for the remote terminal printing and 
160 points for the offline printer. The difference is due to 
the different number of lines per inch. The accuracy of the 
position of a point depends on the range in both the X and Y 
direction. The maximum error in the Y direction is the range 


of Y divided by 200. The maximum error in the X direction is 
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the range of X divided by either 130 for the remote terminal or 
160 for the offline printer. If the X and Y ordinates of two 
Or more points fall on the same nodal point on the graph, only 
one point will be plotted. If the points are from different 
matrices, the matrix hierarchy (see page 34) determines which 


will be plotted. 


D INTEGRATION 

Integration by Simpson's Rule is restricted to equally 
Spaced points in any one direction. The error is of the order 
of H2 for integration in one dimension and of the order of He 
times HS for two dimensions. H, does not have to equal Hy tits 
there has to be an odd number of points in any one direction. 
Integration by Gauss Quadrature for six Gauss ordinates is 


exact for all polynomials of degree eleven or less in one dimen- 


Sion and for CUPOL in two dimensions [Ref. 11] 


Loe POLFIT AND CONTUR 

The polynomial, CUPOL, will pass through the values of 
the function at the nodal points of the grid but may not be a 
true representation in the space between nodal points. If the 
range of values of the function is less than Momo. a contour 


map will not be printed. 
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V. HOW TO CALL THE MATRIX INTERERE PIV oY Solem 


A. TIME SHARING (IMIS VERSION) 

At the Naval Postgraduate School, IMIS is stored in a 
library (AED) with the name, LIT, so that any number of uSers 
can use the system at any one time. In order to call LIT, it 
must be called as a subroutine. A Fortran program called IMIS 
has to be written as follows: 

CALL LIT 

SLOP 

END 

There are two other libraries which must be loaded along 
with AED library in order to use LIT. This can be done with an 
executive program called SMIS as follows: 

FORTRAN IMIS 

GLOBAL LOADER TXTLIB AEDLIB SYSLIB 

LOAD IMIS 

START 

SMIS is then executed with the command, $ SMIS. The 
computer will then type out the executive program with each 
command preceded by a time. Next it will type 

EeoeCULTON BEGINS ... 

Be this time you are in IMIS and ready to begin the first 


problem. The first command will be START. 


i. BATCH PROCESSING (DSMIS VERSION) 
The first card must be a job card for the batch processing 


System as follows: 
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//aaaznnmny,JOB(nnnn, SsssFP,TTTT) , ‘IDENTIFICATION ' 
where 

aaa = The first three letters of the user's last name. 

iz, = Job sequence number 

nnnn= User's number. 

ssss= Project number. 

TTTT= Student section or faculty number. 

The next six cards are control cards. Card five is used 
only when using the punch command. The cards are as follows 
stating in column one: 


//TJOBLIB, DD DSNAME=DSMIS.F59CI , UNIT=2314,DISP=(OLD,KEEP), 


b 


//pbbbbbbbbbbbbY OLUME=SER=MARY 


//ppEXECPPGM=DSMIS , REGION=160K 
//FT06F001,DD,SYSOUT=A 
//¥T07F001,DD,SYSOUT=B 


//¥TO05F001pDDp* 


The a 1s where blanks are left. 


The problem deck follows with the first command a START 
command and the last a STOP command. The last card of the 
deck is a standard delimiter card: 

/* 

On the job request card, check the item marked DISK PACK 

under the SPECIAL INPUT heading and write in "MARY". If the 


punch command is being used, under the heading SPECIAL OUTPUT 


check the item marked PUNCHED CARDS, 
4] v 





Vi. SAMPLE PROBLEMS 


A. THE ELGEN SY¥o2en 


The matrix equation 


[A}) {x} = A[B) {x} 


where 
0.0 -1.0 0.0 
[A] = -1.0 220 -1.0 
0.0 -1.0 Ze0 
and 
1S 0.0 0.0 
fe) = 0.0 2.0 0.0 
0.0 0.0 Ze 


is to be solved. The entire sequence of events is as shown 
below. Everything written in lower case was typed by the 


operator. The computer response iS in upper case. 


&® 
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dOLS 
dois 
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B. CONTOUR MAP 
A contour map for a surface described by the ordinates at 


the sixteen nodal points of the grid shown on page 47 is 


desired. 
10.0 LZ 8.0 10.0 
i) = ZO lero 9.0 8.0 
8.0 9.0 ee 0 12 0 
10.0 8.0 12.0 LOO 
and 
H,=Hy=1 


The results are as follows: 
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The command 
Seftline punch file £f£t07f00l 
will cause the matrix B to be punched onto cards with the 
indicated format of four numbers per card with a field width 
Srel> and 8 Significant figures. 

If the input for the contour map had been ene coefficients 


of CUPOL, the first page of output would have been as follows: 


THIS 1S THE FIRST PAGE OUTPUT COEFFICIENTS AS INPUT. 


=0.333300 X3Y3 1.500000 X3Y2 -2.833000 X3Y 2.000000 X3 
1.500000 X2Y3 ~6.750000 X2Y2 12.750000 X2Y -9.000000 X2 
=2.833000 XY3 12.750000 XY2 -18.750000 XY 9.000000 X 
2.000000 Y3 -9.000000 Y2 9.000000 Y 10.000000 
Q= 7.083 7699378 RANGE ( 0. 7083770 Ol 0.1295°94D 02 ) 
A= 72.3776 B= 7.6713 C= 7.9651 D= 8.2589 E= eso2! 
F= 8.8465 G= 9.1402 H= 9.4340 1= 9.7278 J= 10.0216 
Q= 10.3154 R= 10.60°2 $= 10.9029 T= 11.1967 U= 11.4905 
V= 11.7843 W= 12.0781 X= 12.3719 Y= 12.6656 z= 12.9594 
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oF GRAPH 


A plot of the seune ion 


x = e t9t sin(54.99t) 
is desired from t = 0.0 to t = 0.3. The values were calculated 
in a separate program and the values of X were punched onto 
Cards. They were prepared for offline reading into file 
FTO4FO001. 
The first card is a standard CP67 offline read card as 
foliows starting in column one. 


CP6/7USERID)) GUxxSSSS)),,,p)TTITT 


where 
XxX = Terminal number user will be uSing. 
ssss= User's account number. 
TTTT= User's name. 
The second card is as follows: 


OFFLINE) READ, FILE,FT04F001 


The next card is a LOAD command card with al. in Sl. 
LOAD A 80 2 dag 

The matrix 1S punched onto cards row-wise with the same 
format as one of the options under IMIS and these cards follow 
the LOAD card. AS many matrices as desired may be loaded in 
this manner. The last card is a continue card as follows: 
CONT le 

This card indicates that all matrices have been read 


from file FTO4F001. The deck is then loaded into the computer. 
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After logging onto the time sharing system, the command to read 
the file is as follows: 
offline read * 

OFFLINE READ READ FILE FTO4FO001L 


IMIS is then executed and the sequence is as follows: 
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The command 
OFFLINE PRINTCC FILE FTO8FOOL 


will cause the matrix A to be printed on the offline printer. 
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VII. IMPROVEMENTS 

This system is an open ended system, i.e., any matrix 
operation with fewer than four matrix arguments, four integer 
arguments and two scalar arguments can be coded as a new 
command and added to the system. As coded this program is 
sequential. Each command is executed once after it is read, 
no facilities for automatic repetitive execution of one or more 
commands is included. Such a feature would be very useful if 
1t could be implemented without introducing extraneous 


complexities to the language. 
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READ MATRIX FROM FILE FTO4FOOL 
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